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Abstract. A quantum world-line Monte Carlo method for high-symmetrical quan- 
tum models is proposed. Firstly, based on a representation of a partition function 
using the Matsubara formula, the principle of quantum world-line Monte Carlo 
methods is briefly outlined and a new algorithm using non-binary loops is given for 
quantum models with high symmetry as SU(iV). The algorithm is called non-binary 
loop algorithm because of non-binary loop updatings. Secondary, one example of our 
numerical studies using the non-binary loop updating is shown. It is the problem 
of the ground state of two-dimensional SU(iV) anti-ferromagnets. Our numerical 
study confirms that the ground state in the small N(< 4) case is a magnetic or- 
dered Neel state, but the one in the large N(> 5) case has no magnetic order, and 
it becomes a dimer state. 



1 Introduction 

The quantum world-line Monte Carlo (QMC) method has a long history. It 
has been used in many numerical studies of condensed-matter physics since 
it first proposed in 1970's 1 . However QMC simulations suffered as the 
result of some defects in conventional QMC algorithms. For example, near a 
critical point, or at low temperatures, the correlation time in samples became 
extremely long, thus we could not calculate the canonical ensemble average 
with accuracy. There were same problems in the Monte Carlo simulation of 
a classical system. Fortunately, a new algorithm was proposed by Swendscn 
and Wang, which solves such problems in classical systems. It is called cluster 
algorithm. The idea of the cluster algorithm can apply to quantum cases. 
Indeed, Evertz, Lana and Marcu first proposed a new QMC algorithm based 
on the cluster algorithm, which called loop algorithm [2] . The loop algorithm 
evolves into recent powerful QMC algorithms (See reviews [3] and 0]). 

In the present article, we will focus to the loop algorithm and the gener- 
alization for high-symmetrical quantum models. In Sect. [21 a useful represen- 
tations of a partition function for recent QMC algorithms will be derived. In 
Sect. El the principle of QMC method and the detail of the loop algorithm 
will be reviewed. In Sect. 01 our generalization of the loop algorithm to a non- 
binary loop updating will be proposed. In Sect. [51 we will show one example 
of our numerical studies with the non-binary loop updating: the problem of 
the ground state of the two-dimensional SU(iV) quantum anti-ferromagnets. 
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2 World-line Representations Based on the Matsubara 
Formula 

A density operator has all informations of a quantum system in an environ- 
ment. In particular, if the system is in the heat-bath whose inverse tempera- 
ture is (3, the density operator p(@) is the exponential operator of the system 
Hamiltonian as p(/3) = exp(— where 7i is the Hamiltonian of a quantum 
system. 

The Hamiltonian usually consists of two parts Tio and V. Here TLq and 
V denote diagonal and non-diagonal operators, respectively. The V generally 
represents a perturbation part of the Hamiltonian. Then the non-trivial part 
of the density operator, p'(f3), satisfies a special Bloch equation in which there 
is explicitly only the V operator: 



dp'(P) 
d/3 
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where V(t) denotes the interaction picture of V operator as V(t) = e tn °Ve tn ° 
The solution of |QJ formally satisfies an integral equation as 

p'(0) = I- f ' dtV(t)p'(t), (2) 



Jo 

where an integral variable t is called imaginary time, because the Bloch equa- 
tion is related to the Schrodinger equation through analytic continuation. 

The solution of (0 is written as a series of multiple integrals of the product 
of interaction pictures of V: 

p'{(3) = I - f dhVih) + f dt 2 [ 2 dtj.V{t2)V(ti) (3) 



This representation of a density operator is known as the Matsubara formula. 
The order of interaction pictures V(t) is in descending order. Since the per- 
turbation part V is usually the sum of local interaction Hamiltonians Vb, the 
density operator is a series of multiple integrals of the product of interaction 
pictures of VJ,: 

oo „ 

pW = e -moj2 £ / dt n ---dt 1 (-irv b jt n )---v bl (t 1 ),(4) 

n=0 &„,...,&! ^> tn >...> tl > 

where V = J2b V b- 

The matrix element of operator AB is equal to the summation of products 
of each matrix elements of A and B. 

{V A \AB\V B )= (*a\A\iP}(iI>\B\$ b ), (5) 

{|i/>)}:basis 
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where the set {IV 1 )} is an orthonormal basis in the quantum state space. 
Therefore, inserting an orthonormal basis {|V')} between two operators of 
integrand in the Matsubara formula, the partition function Z can be repre- 
sented as a multiple integral with respect to three kinds of variables, ipi, hi 
and t f. 

Z = Trp(/J) (6) 

OO „ 

= E E E J WnMihibihiU}), (7) 

n=o^ n ,... ) ^i(V'»+i=^i) K*---M p>t n >...>t 1 >o 

n 

W n ({1>i},{bi},{t<}) =e- 0Ho ^l[(^ +1 \{~V bi (U)dU}\^), (8) 

i=l 

where H {ipi) = (tpi\H \^i)- 

The status of a variable set ({i/ii}, {bi}, {U}) is called world-line configu- 
ration, because it can be represented as a set of world- lines. Figure ^ shows 
a world-line configuration of an s = 1/2 Hcisenberg anti-ferromagnet (HAF) 
model on a four-sites chain. The X axis refers the position of spin sites, from 
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Fig. 1. A graphical representation of a world-line configuration for an s — 1/2 
Heisenberg anti-ferromagnet model on a four-sites chain. Solid and dotted lines 
denote spin states +^ and — |, respectively. A position of a waiving line corresponds 
to the status of a variable set (bi,ti). 



1 to 4, and the Y axis refers the imaginary time, from to /?. In this example, 
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the Hamiltonian consists of local interaction Hamiltonians Vi, V% and V3. The 
Vb is the Hamiltonian of the HAF interaction between sites b and 6+1. Then 
a waving line is drawn at the position corresponding to the status of (bi,U). 
The vertical line denotes the status of ipi between two In the s = 1/2 
case, we need two line types, for example, solid and dotted lines which de- 
note up and down spin states of an s = 1/2 spin, respectively. In this way, 
this figure is in one-to-one correspondence with the world-line configuration. 
Usually, the waving line is called vertex and the position at where a local 
spin configuration changes is called kink. Since we will only show examples 
of quantum spin models in the present article, the status of variables tpi is 
called spin configuration. 

From (0), the function W n can be regarded as the Boltzmann weight 
of a world-line configuration. Therefore, using the Matsubara formula, the 
partition function can be transformed from quantum to classical one. 

3 Quantum World-line Monte Carlo Method 

Using the same technique for the partition function, we can get a similar mul- 
tiple integral representation of a canonical ensemble average of an observable 
A. 

(A) S ™M= £ i( ^/n({i}^,Mt i » (9) 

where A(tpi) = (?/>i|A|V>i). Here the symbol {6+ {ti}) denotes the mul- 

tiple integral for variables {{ipi}, {bi}, {ti})- In other words, we sum up with 
respect to all world- line configurations. 

Since Ylu^} {&;} {t ; }) ^™ = tne function W n /Z in © can be regarded 
as the probability of a world-line configuration. Therefore, if we can sample a 
world- line configuration with the probability W n /Z, the canonical ensemble 
average (A) can be calculated as the average of samples of the function A(i/)i). 
This is the principle of the QMC method. 

Many algorithms have been proposed to make a sample of world-line 
configurations. However the conventional algorithm before the loop algorithm 
suffered from some problems. In general, since the probability distribution 
of world-line configurations is complex, a sampling algorithm of world-line 
configurations is based on a Markov process. Samples in a Markov process 
are usually correlated for some time steps which is called correlation time. 
Unfortunately the correlation time in the conventional sampling algorithm 
becomes very large near a critical point and at low temperatures, because 
the size of an updating unit of a world-line configuration is fixed when the 
correlation length of a quantum system increases. Therefore it is difficult to 
calculate a canonical ensemble average of a quantum system at interesting 
points using the conventional algorithm. 
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3.1 Loop Algorithm 

The loop algorithm is one of recent developed QMC algorithms. It uses a 
geometrical object which called loop and updates a world-line configuration 
with loops globally. For some important quantum models, we can prove that 
the size of loop in the loop algorithm is equal to the correlation length of 
the quantum model. This is the reason why the loop algorithm will avoid 
problems in QMC simulations. In fact, the loop algorithm has some good 
properties as follows. The first good property is that the correlation time in 
samples is very short at a critical point and low temperatures. The second one 
is that various types of world-line configurations can be sampled naturally, 
which called grand canonical sampling. In practice, it is difficult to change the 
value of an order parameter without artificial techniques in the conventional 
algorithm. 

The procedure of a loop algorithm consists of three steps. Figure [5] illus- 
trates these steps for an s = 1/2 HAF model on a four-sites chain. The first 
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Fig. 2. A loop updating for an s — 1/2 Heisenberg anti-ferromagnet model on a 
four-sites chain. Solid and dotted line denote spin states +i and — i, respectively. 
The gray loop in the third sub-figure is only flipped. 



step is to update vertexes. In this step, the number of vertexes is changed and 
those positions are moved, but the spin configuration is fixed. We separately 
update vertexes and the spin configuration. The second step is to decompose 
spin variables into loops. In Fig. [21 they are decomposed into five loops. The 
third step is to flip each loop with a probability 1/2. The flipping of a loop is 
to change the sign of spin variables on the loop. In the quantum s = 1/2 spin 
model, a basis of a quantum state space can be made of the direct product 
of two states, up and down states, on each spin site. Thus a spin variable 
takes only two values as ±1/2 and it can be flipped. Although the gray loop 
in Fig. |3 is only flipped, the spin configuration is drastically changed to the 
right one. Since there are great differences between the initial and the final 
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world-line configurations, the correlation time of a loop algorithm is very 
short. 

Since the sampling distribution defined by the weight of a world- line con- 
figuration is very complex and the partition function is unknown in advance, 
we use the Markov process whose equilibrium distribution of configurations 
is equal to the desired sampling distribution. In order to make the Markov 
process, it is enough to satisfy the detail balance condition for any two con- 
figurations A and B in the Markov process: 

W(A)P(A -> B) = W(B)P(B -> A), (10) 

where W( ) is the weight function of a desired distribution and P(A — > B) 
is the transfer probability from A to B. In other words, the density of ran- 
dom walkers is controlled through their transfer probabilities which balance 
with the desired sampling distribution. In the loop algorithm, two updatings 
for vertexes and spin configurations are alternately done and each transfer 
probability satisfies the detail balance condition. 

Firstly, the updating of vertexes is considered. If a vertex V& is inserted at 
an imaginary time t between tk+i and tk (tk+i > t > tk), the weight of a new 
world-line configuration W n +i is equal to the present one (W n in ijHJ) multi- 
plied by a matrix element of V& with dt: (ipk+i \{— Vb(t)dt}\ipk+i) ■ The ratio 
of two weights of pre- and post-inserted a vertex Vt, between the imaginary 
time t and t + dt is proportional to the infinitesimal dt. Therefore, from (|10[1 . 
the accept probability of the inserting of the vertex is proportional to the 
infinitesimal dt. Such stochastic process is called Poisson process. Scanning a 
world- line configuration with the Poisson process, we can insert vertexes un- 
der the detail balance condition. The thing to do next is to remove a vertex. If 
a vertex is not on a kink, the ratio of weights between pre- and post-removed 
configurations is infinite. From (|10|) . we can always remove a vertex on a non- 
kink. On the other hand, if a vertex is on a kink, the ratio of them is exactly 
zero. Thus it can not be can removed. In summary, the stochastic inserting 
process of vertexes is a Poisson process whose intensity is equal to a matrix 
element of a vertex. And the removing a vertex must be done if and only if 
the vertex is not on a kink. 

This updating procedure of vertexes by a Poisson process is applicable 
to various quantum models, and it is employed in various QMC algorithms. 
However, on the contrary, it is difficult to find an efficient updating of a spin 
configuration because of the interaction of spin variables thorough vertexes. 
In particular, the weight of random generated spin configurations is usually 
equal to zero, because the matrix element of each vertex often takes a value 
0. Thus we need an ingenious updating procedure. 

Fortunately, it has been found for some important quantum models. In 
such cases, the interaction Hamiltonian takes the special form which is called 
delta operator. The delta operator is defined as the matrix element taking 
only a value or 1. For example, the HAF interaction Hamiltonian is equal 
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to a special delta operator A(gn) divided by 2 with a special unitary trans- 
formation: 



- H* AF = -U- 1 (S< .Sj-\)U= l -A{g H ). (11) 

Furthermore, the condition of the matrix element taking a value 1 is as fol- 
lows: 

^^1^)1^^ = 10 (otherwise) ' ' (12) 

where <Ji(<Jj) denotes a spin state on a site and it takes an eigenvalue of 
an s = 1/2 z-direction spin operator as <r, = ±4. 

The condition of taking a value 1 for a delta operator can be often rep- 
resented by a graph. If a delta operator of an interaction Hamiltonian has 
the graphical representation, we can construct the loop algorithm on the 
quantum model. For example, the condition of the A(gn) operator is simply 
represented by a horizontal graph in Fig. [31 (a) . The horizontal dotted line 
in Fig. |21 (a) denotes the necessary condition: the sum of two connected spin 
variables must be equal to zero. From the definition, a matrix element of a 
delta operator takes only a value or 1. Therefore the value of a vertex weight 
is not changed, as long as a modification of a spin configuration matches the 
graphical representation. In the horizontal graph, the modification of the two 
connected spin variables is only restricted, thus disconnected spin variables 
are not related to each other: o~i and a[. Furthermore, the flipping of two 
connected spin variables together is always allowed, because it satisfies the 
condition and the weight is not changed. Therefore, we can always flip each 
set of connected spin variables independently. 

In summary, the procedure of updating a spin configuration consists of 
three steps. In the first step each vertex is replaced with a corresponding 
graph. The second step is to decompose spin variables into a set which defined 
by graphs. Since a set of connected spin variables always becomes a loop for 
the s = 1/2 HAF model, this updating algorithm is called loop algorithm. 
The third step is to flip each loop or set randomly with a probability \. The 
loop takes one of two states, because the spin variable takes only a value ±2 
for the s — 1/2 HAF model. Thus we call it binary loop. 

The loop updating was very successful, because the size of loop is related 
to the physical correlation length of a quantum model. And by the loop algo- 
rithm, many important results were numerically found for various quantum 
models. 

Although we showed only a quantum s = 1/2 spin model, the loop 
updating is also applicable to quantum large spin models. For such cases, 
Kawashima and Gubernatis proposed a mapping to a quantum s = 1/2 spin 
model which called split-spin technique jSJ. The idea is that an s — m spin 
operator is replaced by the sum of 2m s = 1/2 spin operators with a per- 
mutation operator. Since they are isomorphic, we can simulate a quantum 
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s = m model as a multi quantum s = 1/2 model with binary loop updatings. 
However, since the size of the new quantum space becomes large, we need 
more memories and complex treatments. In the next section, we will propose 
a new loop updating which acts on an original quantum space directly. 

4 Non-binary Loop Algorithm 

The symmetry of a model is related to the graphical representation of a 
vertex. For example, the HAF model has the SU(2) symmetry and the HAF 
interaction term is represented as the binary horizontal graph. Fortunately, 
there are such special graphical representations in high-symmetrical cases. 

The s = 1 bi-linear bi-quadratic (BLBQ) model is interesting theoret- 
ically and experimentally [S], because it has some integrable points for a 
one-dimensional case and it is related to the effective model of Na atoms in 
an optical lattice [7|. The Hamiltonian is as follows: 

Hf^ BQ = (cos 6)Si ■ Sj + (sm6) (S, ■ Sj) 2 . (13) 

The BLBQ model has SU(3) symmetry at special points, for example, = — ? 
and — The BLBQ interaction terms at SU(3) points have special graphical 
representations whose states are non-binary as and ±1. For example, the 
BLBQ interaction term at the SU(3) point = — f becomes a delta operator 
and the condition taking a value 1 is same to that in an s — 1/2 HAF case, 
except for a spin variable takes three values, ±1 and 0. As a result, we can 
use the same horizontal graph in Fig. [31 (a) for the BLBQ interaction term. 
However, the loop, which defined by the SU(3) horizontal graph, takes three 
states, not binary states. 

We can generalize the horizontal graph in the SU(iV) case. In the SU(A r ) 
case, since the spin variable takes N values (ct, = , • • • , ( ' Af ~ 1 ^ ), the 

number of loop states which defined by the SU(iV) horizontal graph is N. 
Thus we call it non-binary loop. Another type of graph can also be generalized 
in the SU(AT). For example, we found the generalized cross graph gc in Fig-El 
(b) which corresponds to the other SU(3) point 9 = — 

wi^>iw - {J j^f"* -''-•». <») 

It is easy to check that generalized graphs have the SU(iV) symmetry. 

The loop defined by the SU(iV) graph has N possible states of equal 
weight, because the exchange of any two states corresponds to one of SU(iV") 
transformations. Therefore we can choose one of N possible states for each 
loop with equal probability. 

In summary, the procedure of the loop updating using the SXJ(N) graph 
is similar to that of the conventional loop algorithm with the binary graph. 
In the first step each vertex is replaced by a corresponding SU(A r ) graph. 
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And the second step is to decompose spin variables into a set which defined 
by the SU(iV) graphs. In the third step we choose one of N possible states for 
each set with equal probability. The only difference between the conventional 
and new loop algorithms is the number of possible loop states. 




Cross graph 



Fig. 3. Graphs for various types of vertexes and a non-binary loop updating for an 
s — 1 BLBQ model on a four-sites chain, (a) Horizontal graph, (b) Cross graph, 
(c) The non-binary loop updating for an s = 1 BLBQ model at 9 — — ~. Solid, 
dashed, and dotted lines denote the spin states, +1, 0, and -1, respectively. The 
gray loop in the second sub-figure is only changed. 

The loop algorithm with non-binary loops is called non-binary loop al- 
gorithm @|. It is very useful to study high-symmetrical quantum models, 
because of the simplicity. Without a split-spin technique, we can directly 
simulate the SU(iV) model in the original quantum space. The relation be- 
tween the conventional and new loop algorithms is similar to that between 
Swendsen-Wang algorithms for Ising and Potts models. 

5 Application: Ground States of Two-Dimensional 
SU(iV) Quantum Anti-ferromagnets 

Some quantum models already have been simulated with the non-binary loop 
algorithm. In the present section, we will focus on the numerical study of 
ground states of two-dimensional SU(iV) quantum anti-ferromagnets. 

The existence of a short-range resonating valence bond (RVB) spin-liquid 
is a very interesting problem for low-dimensional quantum systems. A RVB 
spin-liquid exhibits a finite gap for spin-excitations and it has only short- 
range order and it does not break any lattice symmetry. Although the RVB 
spin-liquid may be related to the mechanism of the copper-oxide supercon- 
ductors, the ground state of the two-dimensional s = 1/2 HAF model, which 
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is an effective quantum spin model of copper-oxide superconductor, is not 
a RVB state, but it has a magnetic order which called Neel order. As the 
spin s is decreased, the quantum fluctuation decreases the Neel order, but 
the ground state remains ordered even for s = 1/2. Another route to increase 
quantum fluctuations is in the symmetry of quantum models. We studied 
the SU(A) quantum anti-ferromagnets which is one of generalizations of the 
HAF model and which has the higher symmetry than SU(2). 

Our Hamiltonian of SU(iV) quantum anti-ferromagnets is defined by the 
generator of SU(AT) algebra as the HAF model: 



a, (3=1 

[J$W, ■#(?)] = % (%Jg(i) ~ %J?(ij) > (16) 

where Jg(i) is an SU(AT) generator on site i. In the SU(2) case, this Hamilto- 
nian is equal to the HAF model. In our anti-ferromagnetic SU(iV) model, the 
representation of SU(iV) generator is a A-dimensional matrix which called 
fundamental representation and the matrix on one sub-lattice is conjugate 
to that on the other sub-lattice. In this representation, the model can be 
expressed in terms of SU(2) spins with s = (N — l)/2. If the SU(AT) Hamil- 
tonian applies to the \nm) state, the result equals the superposition of all \ll) 
states divided by N: 

H™ [N) \nm) = -8 nrn [^ £ = -S nrn \(ij)). (l = -l), (17) 

V l=-a, — ,a J 

where \nm) denotes the simultaneous eigenstatc of z-direction SU(2) spin 
operators on sites % and j. For convenience, the special superposition state 
will be called (ij)-bond state in the following. 
The (12)-bond state is a ground state of a two-sites case. And it has no 
magnetic order, because it is the superposition of different magnetic ordered 
states. In a four-sites bond state | ( 12) (34)) is an eigenstate of Tif^ 

and W34 whose eigenvalue is -1. However, if the interaction Hamiltonian 
does not match the covering of dimers over sites in a bond state, the state 
converges to a zero-eigenstate when N becomes infinity: 

^ 3 U(JV) |(12)(34)) = -1|(41)(23))^0 (A^co). (18) 

Therefore, in the large N limit, the ground state maximizes the number 
of nearest-neighbor bonds, because the eigenvalue of a whole Hamiltonian is 
minus of the number of nearest-neighbor bonds. In short, in the large N limit, 
due to strong quantum fluctuations, the ground state has no magnetic order. 
However, in a one-dimensional chain, the number of maximum bond state 
is only two Since they break the translational symmetry, they are not a 
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RVB state, which called dimer state. Read and Sachdev studied the ground 
state of SU(iV) models on a two-dimensional square lattice theoretically 8 . 
Their conclusion is also that the ground state on a two-dimensional square 
lattice is a dimer state in the large N limits. 

In the N — 2 case, the ground state is a magnetic ordered Neel state. On 
the other hand, in the large ./V limit, it has no magnetic order and it becomes 
a dimer state by Read and Sachdev's prediction. However the intermediate 
./V cases remain problems. The numerical study by Santoro et al. suggested a 
RVB spin- liquid ground state in the A" = 4 case jEJ]- The model of Santoro 
et al. is equal to our SU(4) model and it also appears as a special point in 
coupled spin-orbital models. Santoro et al. used the Green function Monte 
Carlo method in their numerical simulations. Due to the computational times 
by their numerical method, the size of their simulated lattices was limited to 
only L = 12. 

Since the SU(A') interaction Hamiltonian is equal to the delta operator 
A(gn) represented by the non-binary horizontal graph, this model can be 
simulated for various N and large square lattices by a non-binary loop al- 
gorithm. Indeed, we explored the system size L up to 128 and N up to 8 
Our simulations have been performed at low enough temperatures to 
be effectively in the ground state. The lowest temperature for N — 4 and 
L = 128, for example, is j^g. 

In order to search the order of a ground state, we have measured two 
quantities which are related to magnetic and dimer orders, respectively. The 
one is a static structure factor S(ir, n) of staggered magnetization m s and 
the other is a dimer static structure factor S b (tt,0): 



%^ = ((m s ) 2 ) (m s ee -1 $>irS? I (19) 

and 

= ((<W)) 2 > (^o) = £(-l) r " W+J ■ (20) 

If an order exists, the related structure factor divided by L 2 converges to a 
finite value in the large- L limits, but if the order is absent, it decreases to 
zero as 1/L 2 . 

Figure |U shows a clear evidence for a Neel order in JV < 4, because the 
spin structure factor divided by L 2 converges to a finite value for N < 4. 
Therefore the ground state is not a RVB spin-liquid even in the N = 4 
case. On the contrary, in the N > 5 cases, the spin structure factor divided 
by L 2 decreases as 1/L 2 . Therefore there is no Neel order in JV > 5. The 
disagreement of our conclusion for the A^ = 4 case with that of Santoro et al. 
is not in the raw numerical data. Their data are only in the L < 12 region. 
Therefore, the disagreement is solely due to the small system sizes studied in 
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Fig. 4. Static structure factors Sr^*) for 2 < TV < 8. The straight line representing 
the power law L~ 2 is drawn for comparison. Estimated statistical errors are not 
shown, because they are equal to or smaller than the symbol size. The inset presents 
the data for TV = 3, 4, 5, 6 in the linear scale, together with the best fitting curves 
obtained by the method of least squares. (Adopted from Harada and Kawashima 

ED) 



On the other hand, Fig. [5] shows a clear evidence for a dimer order in 
N > 5. While the dimer structure factor divided by L 2 decreases to zero in 
N = 4, it converges to a finite value in N > 5. Therefore the dimer order 
exists in N > 5 and the ground state is not a RVB spin-liquid state. It is 
consistent with the theoretical prediction by Read and Sachdev in the large 
N limits. 

In summary, from QMC simulations for large square lattices by the loop 
algorithm with non-binary loops, we confirmed a direct transition from the 
Neel ground state (N < 4) to the dimer ground state (N > 5) for two- 
dimensional SU(iV) quantum anti-ferromagnets. These comprehensive nu- 
merical studies can be effectively done by the loop algorithm with non-binary 
loops. 



6 Conclusion 

We have introduced the QMC method with non-binary updatings for quan- 
tum models with high symmetry. For example, using graphical representa- 
tions of the interaction Hamiltonian, we can easily simulate a high-symmetrical 




1/L 

Fig. 5. The k = (ir, 0) dimer structure factors ) for N = 4, 5, 6 in logarithmic 
scale. The inset is the linear-scale plot. The solid lines in the inset are the best 
fitting curves of least squares based on the L > 8 data. (Adopted from Harada and 
Kawashima |llp 

quantum model by a loop algorithm with non-binary loops. We have only 
focused to the loop algorithm, but the idea of a non-binary updating for 
quantum models with high symmetry can be combined with the other algo- 
rithms as the worm or directed loop. We have shown one example for using 
a non-binary loop algorithm. It is the problem of the ground state of two- 
dimensional SU(iV) quantum anti-ferromagnets. Our conclusion is that our 
SU(N) models have no RVB spin-liquid ground state. 
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